Analytical expressions for primary Bjerknes force on inertial cavitation bubbles. 
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The primary Bjerknes force is responsible for the quick translational motion of radially oscillating 
bubbles in a sound field. The problem is classical in the case of small-amplitude oscillations, for 
which an analytical expression of the force can be easily obtained, and predicts attraction of sub- 
resonant bubbles by pressure antinodes. But for high-amplitude sound fields, the bubbles undergo 
large amplitude nonlinear oscillations, so that no analytical expression of the force is available in 
this case. The b ubble dynamics is approximated on physical grounds, following the method of 
iHilgenfeldt et all [J. Fluid Mech., 365, 171 (1998)], but carefully accounting for surface tension. 
The analytical expression of the maximum radius of the bubble is recovered, the time of maximum 
expansion is noticeably refined, and an estimation of the collapse-time is found. An analytical 
expression for the time-varying bubble volume is deduced, and the Bjerknes force is obtained in 
closed form. The result is valid for any shape of the sound field, including purely standing or purely 
traveling waves, and is ready to use in a theoretical model of bubble clouds evolution. Besides, the 
well-known sign inversion of the Bjerknes force for large standing waves is recovered and the inversion 
threshold in the parameter space is obtained analytically. The results are in good agreement with 
numerical simulation and allow a quantitative assessment of the physical parameters effect. It 
is found that either reducing surface tension, or increasing the static pressure, should produce a 
widening of the bubble-free region near high-amplitude pressure antinodes. 

PACS numbers: 47.55.dd, 43.35.Ei 



I. INTRODUCTION 



When excited by a sinusoidal sound field, gas bubbles 
undergo radial oscillations. Most of the practical appli- 
cations of this phenomenon, known as acoustic cavita- 
tion, use high-amplitude sound fields, of typical ampli- 
tude greater than the static pressure, so that the liq- 
uid is under tension for some part of the cycle. In such 
conditions, whatever the frequency, two distinct dynamic 
bubble behaviors can be clearly divided by the so-called 
Blake threshold [ll-ll]' in the tension phase, very small 
bubbles are retained to grow by surface tension. Con- 
versely, larger ones suffer an explosive expansion followed 
by a violent collapse, responsible for chemical [Sj, me- 
chanical effects d, 0| and sonoluminescence [8l-[io|. The 
latter oscillation regime is known as "inertial cavitation" . 

Bubbles in liquids experience various hydrodynamic 
forces. The buoyancy force is the most familiar one, and 
is the pressure force that an sphere of liquid replacing 
the bubble would experience. This remains true in a ac- 
celerating liquid ll[ , and the generalized buoyancy force 
experienced by the bubble is —WV where V{r, t) is the 
pressure that would exist at the center of the bubble if it 
were absent, and V{t) the bubble volume. For a bubble 
oscillating radially in a sound field, both Vir^t) and V{t) 
are oscillatory quantities so that the time-average of the 
product over one cycle is not zero. The bubbles experi- 
ences therefore a net force known as "primary Bjerknes 



force" hahi 
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(1) 
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The Bjerknes force can be easily calculated from the 
knowledge of both the shape of the sound field and the 
bubble dynamics. A classical result is that for low- 
amplitude standing waves, sub-resonant bubbles are at- 
tracted by pressure antinodeSj_while bubbles larger than 
resonant size are repelled |14l4l6| . For the case of strong 
driving pressures, sub-resonant inertial bubbles can also 
be attracted by pressure antinodes, which constitutes the 
basic principle of SBSL levitation cells How- 
ever it has been shown by numerical calculations that 
above a given threshold, the primary Bjerknes force on 
sub-resonant inertial bubbles undergoes a sign change 
[l8| . This behavior is due to the resonance- like re- 
sponse curve (termed as "giant resonance" by Lauterborn 
and co-workers [19]) of the bubble just above the Blake 
threshold, which is a physical consequence of the effect 
of surface tension. Experiments indeed demonstrate that 
above a certain driving level, no bubbles are visible in the 
neighborhood of large pressure antinodes [l^ . 

Quantitative agreement between theory and experi- 
ment has been found in the case of linear or q uasi-linear 
oscillations [1^1 ■ Particle simulations [13, l2l| were also 
found in excellent agreement with recent experiments in- 
volving inertial bubbles [13|- While the Bjerknes force 
can be calculated analytically for linear bubble oscilla- 
tions, only numerical results can yet be found for inertial 
bubbles [3 [13 • analytical expression for the latter 
would first be helpful in particle or continuum models, 
describing the self-organization of bubbles, in order to 
get more efficient calculations. Furthermore, analytical 
results allow a direct assessment of the sensitivity of the 
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force to the physical parameters, and the estabhshment of 
scahng laws. These two objectives raotivated this study. 

Owing to the strong nonlinearity of the bubble dynam- 
ics equations, inertial cavitation has long been thought 
intractable analytically, up to the seminal papers of 
Lofstedt et al. 24 1 and Hilgenfeldt et al. who demon- 
strated that several terms of the Rayleigh-Plesset equa- 
tion (RP) could be neglected during the explosive ex- 
pansion of the bubble. This theoretical breakthrough 
allowed to obtain scaling laws for the maximum radius 
of the bubble and the time of maximum expansion. In 
this paper, we closely follow the approach of Hilgenfeldt 
et al. [l| and refine their analytical solutions in order to 
account more precisely for the effect of surface tension. 
The approximate dynamics found are then used to ob- 
tain an analytical expression of the bubble volume. The 
latter are then conveniently recast in order to obtain the 
Bjerknes force ([1]) in closed form, in any acoustic field, 
including the two extreme cases of traveling and standing 
waves. Finally, in the latter case, we seek an approximate 
expression of the Bjerknes force inversion threshold, evi- 
dencing the role of surface tension. 



Rayleigh-Plesset (RP) equation [ll,!!!-!!^: 
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where po is the hydrostatic pressure, Pg {t) is the gas pres- 
sure, p, and c/ are the density, viscosity and sound 
speed of the liquid, respectively, and a is the surface ten- 
sion. The ambient radius of the bubble Rq is the radius 
that would have the gas the in absence of the sound field. 

Time is non-dimensionalized by the angular frequency 
w, and in order to obtain a formulation consistent with 
Ref. [H, we set 



Pa + 'P{r,t) =po{l -pcosx). 



so that 
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Ult + <l){r) ~~ TT. 
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II. PRIMARY BJERKNES FORCE 
A. Acoustic field 

We assume that the acoustic field in the liquid is mono- 
harmonic at angular frequency w, and defined in any 
point r by 



-p(r, t) = Pir) cos [ujt + (j){r)] 



(2) 



This expression may represent a traveling wave, a stand- 
ing wave, or any combination of both. We also define the 
pressure gradient in general form as 



dV 

— (r,i) = Giir) cos[ujt + %l:i{r)] 



(3) 



Using X as the time- variable, and non-dimensionalizing 
pressure with poj equation (j4|) can be written as: 



RR" + £i?'2 
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p cos X — 1 



where primed variables denotes d/dx, 

Rros = u;-\3po/pY/^ 
is the resonance radius, and 

as = 2a/poRo 



(8) 



(9) 



(10) 



where the fields Gi and tpi can be expressed as functions 
of P and (f) once the acoustic field is known. The following 
two extreme cases deserve special consideration: 

• for a standing wave, (/)(r) = (po, so that Gi{r) — 
dP/dxi and ^/'i(r) = 0o, 

• for a traveling wave, P{r) ~ Po and 0(r) = —k.r 
so that Gi{r) = fc^Po and ipi{r) = ip{r) — tt/2. 



B. Bubble model 

The radial oscillations of a gas bubble in a liquid under 
the action of the sound field can be described by the 



is the dimensionless Laplace tension. 

Several models can be used for the bubble internal 
pressure Pg, [28l - [3^ . As will be seen below, we are 
mainly interested here in the expansion phase of the bub- 
ble, during which the density of the gas in the bubble 
remains weak, so that the precise choice of the thermal 
bubble interior's model is unimportant. However, in or- 
der to assess the validity of the approximate expressions 
developed hereafter, simulations will be performed by us- 
ing the Keller equation [s^, [s^]. The bubble interior 
is modeled by using a thermal diffusion layer following 
Ref. [13, neglecting water evaporation and condensation 
through the bubble wall. In the remaining part of the 
paper, we will consider air bubbles in water (cr — 0.072 
N.m^i, Po = 101300 Pa, p = 1000 kg.m^^ q = 1498 
m.s^^, fj, = 10^^ Pa.s). 
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C. The Bjerknes force 

The primary Bjerknes force acting on a bubble is de- 
fined as 



{V{t)VV) , 



(11) 



wliere V{t) is the instantaneous bubble volume. The av- 
erage is taken over one acoustic period, so that, using 
Eq. 0: 



-G^{r)- I V{t)co8[Lot + ^,{r)] dt. (12) 



Using the dimensionless time x defined by ([7]) and the 
periodicity of V , the latter expression becomes 

1 /-^^ 

FBi = Gi{r)— / V {x) cos.[x - (j){r)+ipi{r)] dx. 
Stt Jo 

(13) 

The generic problem is therefore to obtain an approx- 
imate analytical expression for the integral 
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V {x) cos (x — xa) dx, 



(14) 



valid for any bubble dynamics, and for any value of xq. 
The problem can be easily solved for small-amplitude lin- 
ear oscillations [iBl- Here, we focus on the case of inertial 
oscillations, that is for any combination of parameters 
(p, i?o) above the Blake threshold. The special cases of 
standing waves and traveling waves can be simply recov- 
ered by setting, respectively, xq = and xq = 7r/2. 



III. APPROXIMATE EXPRESSIONS 

A. Bubble radius 

The method used to obtain analytical formula for the 
bubble radius are mainly inspired from the approach of 
Hilgenfeldt et al. [H. For self-consistency, we will recall 
in this section the main lines of the method, and, where 
convenient, specify the refinements obtained by our ap- 
proach. 

Figure [T] displays the dimensionless bubble radius ([TJi), 
bubble volume (Hb), and driving pressure ([TJ;) in a typical 
case of inertial cavitation of (/ — 20 kHz, i?o = 3 /im, 
and p — 1.4). With the choice of the dimensionless time- 
variable Eq. ([5]), a; = represents the time of maximum 
tension of the liquid. We set 



1 

a;+ — acos — , 
P 



(15) 



and we denote by Xm the time of maximum expansion of 
the bubble, and by Xc the time of its maximum compres- 
sion (see Fig. [TJ. 

It is shown in Ref. [l| that, during the expansion phase 
and most of the collapse phase, the dominant terms in 
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FIG. 1: (a) Dimensionless bubble radius R/Ro; (b) Dimen- 
sionless bubble volume (i?/i?o)^; (c) Dimensionless driving 
pressure 1 — pcosa;. The case considered is a 3 fim air bubble 
in water and p — 1.4. The times —x+ and x+ are the two in- 
stants of zero-crossing of the driving pressure, Xm is the time 
of maximum expansion of the bubble, and Xc the time of max- 
imum compression. The dashed curve in (a) represents the 
approximate dynamics given by Eqs. (|21l) . (|22|) . The dashed 
line in (b) is the final approximation of the bubble volume 
([36l)-(|42|. 



the right-hand-side of Rayleigh equation are the driving 
term p cos x—1 and also the surface tension term asRo/R 
for ambient radii just above the Blake threshold. Fol- 
lowing Ref. [H, we neglect the dependence of the surface 
tension term in R , and replace asRo/R by as/K{p), 
where K{p) will be determined later. The approximate 
Rayleigh equation becomes 



RR" + 



2 

res 



p cos X 



as 
Kip) 



We set, for further use 

A=l 

Besides, noting that 

RR" + : 



as 



K{p)- 
d^R^) 



1/2 



dx"^ 



(16) 



(17) 



(18) 



the right-hand side of Rayleigh equation can be written 
in two different forms: 



RR" 



2^ = - 



l_(P{R^ 1 
2 ^2 
3 d'(i?^) 1 
4 2 



-R' 
RR. 
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Numerical simulations show that B? ^ RR on the 
interval [— a:_i_, while R? ^ RR holds on the inter- 



val 



ff 



Additionally, we found that the latter 



property still holds in fact during almost all the collapse, 
except in its ultimate phase, where the gas and acoustic 
terms become significant again. This could be expected 
since the main part of the collapse is inertially driven 
and that R becomes significant only when the liquid has 
acquired enough kinetic energy. We therefore obtain the 
following equations for the bubble radius, over the inter- 
val [—x+, x^: 



^^2 = g^rcs(PCOSX-^) On , 2;+] , (19) 



<P[R^) 



- R^^g {p COS X — A) on[x+,Xc]- (20) 



These equations are the same as the ones of Hilgenfeldt 
et al. except that the validity of the second is ex- 
tended up to Xc- The first equation can be solved with 
the initial condition R{—x+) — C^OjWhere C ~ 1.6 and 
R{—x+) ~ R{—x+) [1]. The second equation is solved 
by requiring continuity of R{x) and R'{x) at a; = 
Integrating both equations twice, we obtain: 



R^_{x) = -R 



1 — p cos X + p{x + x+) sin x+ 



A 
~ 2" 

+ eRl [l + 2{x + x+)], 



(x + x+) 



(21) 



and 



Rlix) = 



2 

3^rcs 

A 



I — p cos X + p + x^^ sin x+ 



2 ,x^+xl + -x+x 



+ C^Rl [l + 2{x + x+)] 



(22) 



The point (xm,^max) of maximum expansion is ob- 
tained by setting d{R\)/dx = 0, so that Xm is given in 
implicit form by 



psivLXm - Xm + - (j3sinx+ - x+) 



as 
K{p) 



and i?max reads 



^max — RqUPj Xn 



3C 



Rrcs 



(23) 



0, 



Rr- 



2 as 



3K{p) 



h{p,Xm) 



(24) 



where 

f{p,Xra) = C^l+2{Xm+X+)] 
g{P:Xm) = ^ 



(25) 



1 — p cos Xm + P ( — — h x^ J sm x+ 



-2 1^. 



(26) 
(27) 



In order to obtain Xm, the implicit equation (1^ should 
be solved. To avoid this, Hilgenfeldt and co-workers [l| 
developed this equation near 7r/2 at first order, neglecting 
on the one hand as/K{p), and also (i?o/^rcs)^, which 
is appropriate for driving the bubble at low frequencies. 
They obtain 



Cmo =P + ^ ipsmx+ -x+), 



(28) 



which can be further simplified as Xm — P, if P is small 
enough. Plugging the latter into Eqs. ([M|) - (l?7l) . they 
obtain an expression of -Rmax which depends on Rq only 
through the as term in (|24p . The expression of K{p) is 
then determined by using the fact, confirmed numerically, 
that the maximum of the response curve {Rma.^/ Ro){Ro) 
is obtained for an ambient radius Rq very close to the 
Blake threshold 



d ( Rma,x(p,Ro) 



dRo 



Rq 



for Rq = R^ 





4V3 a 



1 



9 po P - 1 



(29) 



This scheme yields a good approximation for i?niax: which 
was the main objective of Hilgenfeldt and co-workers [T], 
but the approximation (|28p of Xm yie-^'is ^ rather large 
error (see dotted line in Fig. [5] and [41|). Since the value 
of the integral is found to be very sensitive to the 
precise location of Xm , we seek a better approximation. 

We therefore revert to the original equations ([^ -(P? ]) . 
The main difficulty lies in the presence of the as term 
in (j23p . which makes rigorously Xm a function of both 
p and Rq. Thus -Rmax depends not only on Rq through 
as but also through x^ in the expressions of /, g and h. 
The condition ([Ml) therefore becomes more complex, and 
should be solved simultaneously with (j23p . We initially 
followed this complex process, but finally found that a 
better approximation of Xm could be obtained by using 
a simple trick. First, as was done in Ref. [s^ we neglect 
the as term in ((23)) and develop the latter near 7r/2, but 
up to second order : 



TT 1 1 , 
2 p p 



TT 



2p 



Rq 
Ric 



1/2 



(30) 
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For low frequency driving, Rq <^ Rrcs, and depends 
only slightly on Rq. We then plug ® in Eqs (Pil) - (P7)) 
and express the condition (|29|) . neglecting dxm/dRo, to 
obtain 



Klip) 



4 



4^/3 



(31) 



We now expand again (1^ near 7r/2, but keeping the as 
term, in which we set K — Ki (p) , to obtain 



TT Ai 



P P 



Al 



2p 



3C 



1/2 



(32) 





where 



= 1 



Q!5 

/^l(rt' 



(33) 



Setting Al = 1 in Xm2^ that is, neglecting the effect of 
surface tension, the result of Ref. |3a, Eq. (pO| . is recov- 
ered. 

Figure [2] represents the variations of Xm for a bub- 
ble of ambient radius Rq — 1 /im (fig. Oa) and Rq — 3 
/ini (fig. [21b) in water. The thick solid lines are the ex- 
act value obtained numerically, and the thin solid lines 
represent Xm2- The agreement is seen to be excellent, 
although a noticeable difference can be seen for Rq = 1 
/Ltm, which originates from the over-simplification done 
when accounting for surface tension by the simple term 
as/K{p) in Eq. (1161) . Also shown is the approximation 
Xmi (dash-dotted line), which does not take surface ten- 
sion into account. This clearly introduces an noticeable 
error on Xm, reasonably corrected by Eq. p2p . Finally, 
the approximation Xm — P proposed in Ref. [l[ is dis- 
played (dotted line). 

Finally, the approximation of i?max can then easily 
be obtained by plugging an a ppr oximation of Xm into 
Eq. p4|) . This was done in Ref. |35| using Xmi , and an ex- 
cellent agreement was found. The gain brought by using 
Xm2 instead of Xmi in ((Ml) remains unimportant, and for 
brevity, we do not present the comparison between the 
analytical and numerical expressions of -Rmax here. 



B. Bubble volume 



FIG. 2: Thick solid line: Xm calculated from numerical so- 
lutions of equation (jU; Thin solid line: a;™, from Eq. (I32|l : 
Dash-dotted line: from Eq. (|30p (ref. [35]); dotted line : 
Xm = p (ref. 1]). Thick dashed line: Xc calculated from nu- 
merical solutions of equation Q; Thin dashed line: Xc from 
Eq. (|43|) . The results are calculated for a bubble of ambient 
radius Ro — 1 fim (a), and Ro = 3 fim (b). 



1. Approximate expression on [—x+,x+] 

Numerical simulations demonstrate that R- is almost 
linear between and x+ (Fig. [Ua), which suggests that 
Eq. (PT|) is almost a perfect square in this interval. We 
then develop the cos term in ([21]) near a; = up to second 
order, and write the result as 

psinx+ — Ax^'' ^ 



' pxj^ sma;^ 



X - 
A 



p- A 
l-p 



(34) 



(psinx+ — Ax+Y 



2{p-A) 

For _R_ to be linear in x, the constant term in the bracket 
must be negligible, so that i?_ can be simplified as 



p - A 



Ax 



+ — p sm x+ 
p - A 



(35) 



The expression of the bubble volume on [—x+, x^] there- 
fore reads 



Approximations of the bubble volume could readily be 
obtained from the approximations ([21]) , ([22]) of the bub- 
ble radius. However, such expressions do not yield ana- 
lytical expressions of the integral ([H)) in closed form, and 
further approximations are therefore required. First, we 
consider frequencies low enough to have R^cs 3> Rq, so 
that the C term can be safely neglected in equations ([^T]) - 



where 



2o .IP- A 1 „ ,3 



(x-xl)^ (36) 



Xl 



Ax+ — p sin x+ 
p- A '■ 



(37) 



which allows to calculate integral ([Ti)) in closed form. 
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2. Approximate expression on [x+,Xc] 



IV. BJERKNES FORCE 



Using equations ([23]) and (|24|) . it can be easily checked 
that, setting y — x — Xm, the expression (1^ of i?+ can 
be recast as: 



p2 d2 



(38) 



where 



i(y) = 2pcosXm sin^ | - 

+psinx„i (sin?/ - y) . 

The bubble volume on , Xc] becomes therefore 

-1 3/2 



L{V) 



(39) 



(40) 



which unfortunately does not yield an explicit integration 
of (fT4|l . Further progress can be done by noting that, 
from Eq. (1^ . i?res and i?max are of the same order of 
magnitude, and that from ([39)) . L{y) — 0{y'^) near y = 0. 
Equation (j40p can therefore be approximated by 



V+ ^ R 



3 

max 



Rr, 



Rn 



L{y) + 0{y' 



(41) 



Thus, to the same order of approximation, L{y) can 
be replaced by any equivalent expression up to order 
4 in y, and the choice must be directed by the ability 
of Vj^ cos{x — Xq) to be integrable in closed form. We 
therefore choose to set = sin^(y/2) + 0{y'^) and 

siny - y = -l/6sin^y + O(y^) in Eq. to finally 

obtain: 



V+{x) 



p p^ 



2{p cos x„i — A) sin^ — 



--psma;„ism y 
6 



O(y^), 



(42) 



which can now yield an explicit expression for integral 

m- 

It can further be noted that neglecting the sin"^ y term 
in the square bracket, and setting sin^ | ~ 2/^/4, V+ is 
found to be zero for 



Rn 



Rr. 



A — p cos x„ 



1/2 



(43) 



A. Analytical expression 

With the expressions of the bubble volume 
and at hand, the integral (|14l) can be calculated 
in analytical form, keeping the contribution of the inte- 
grand only in the intervals [0,a;+] and [a;+,Xc], since V 
can be neglected in the other regions (see Fig. [TJb) . The 
integral is thus the sum of the two contributions: 



where 



I- 



V- (x) cos {x — xq) dx 



V+(a;) cos (x — xo) dx. 



(44) 



Using the approximate expressions (15^ and (|^^ of the 
bubble volume, integration yields 



" p3 



p- A 



3/2 



Ax{Ax'^ — 6) sin(a::+ — xq) 
+3(Ax^ - 2) cos(a;+ - xq) 
+xi{6 — Xi) sinxo + 3(2 — xf) cosxp 



,(45) 



with 

Ax = x^ 
The contribution /+ reads 



Xl. 



1+ — R„ 



[sin(a;c — xq) — sin(a::+ — Xq) 



-Rmaxi?rcs | ^ (P COS - A) 

x[/2(yc)-/2(y+)] 
1 



psmx„, [faiVc) - f3iy+)] >, (46) 



where 



f2{y) 
hiv) 



192 



4sin(y - yo) - sin(2y - yo) - 2ycosyo, 
2 cos(2y + yo) + cos(4y - yo) 
+12ysinyo - 6cos(2y - yo). 



and 



which constitutes a simple approximation of the collapse- 
time. The comparison between this expression and the 
exact instant of minimum radius is visible in Fig. [5] 
(dashed lines). Here again, an excellent agreement is 
found, but deteriorates toward small bubble radii. 



Va^ Xq - x,n, y+^x+~x,, 



The value of / from is displayed in figure [3] 

(thick lines) for i?o = 1 /xm (Fig. [3^), 3 /xm (Fig. [3]d), and 
6 /xm (Fig. ISt), in the case of a standing wave [xq = 0), 
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for drivings ranging from the Blake threshold to p = 2.5. 
In order to get a clear picture, / is drawn in logarithmic 
scale, the solid part of the curves representing a posi- 
tive sign and the dashed part a negative sign. The thin 
lines are the results obtained by solving Q and calcu- 
lating (|14p numerically, for / = 20 kHz. It is seen that 
an excellent agreement is obtained, except for Rq = 1 
/im (Fig [3^) . Particularly, the point of inversion of the 
Bjerknes force is shifted toward large drivings. This fea- 
ture originates from the errors induced on the values of 
Xjn, Xc (see Fig. and R^ax for small ambient radii, 
by replacing the surface tension in the RP equation by 
as /K{p) in (IT6)) . It should be noticed that even the small 
errors visible on the curves of Fig. [Sja yields large differ- 
ences on the estimation of /. This could be expected since 
the phase between V and cos(x — xq) crucially influences 
the value of integral /. 




FIG. 3: Thick lines: value of I predicted by approximation 
(|44|l - (|46| l. for a;o = (standing wave), for an air bubble in 
water, of radius Ro — 1 fim (a), 3 /im (b) and 6 fim (c). The 
solid parts of the curves correspond to / > and the dashed 
parts to / < 0. 



B. Bjerknes force inversion threshold in standing 
waves 

We consider the case of a standing wave a;o = 0, and 
look for an approximate locus in the parameter space 
where the Bjerknes force changes sign. Summing equa- 



tions p5|) and 
for 



where 



46)) . it is seen that integral (|14l) is zero 



(47) 



aaX-^ + aiX + uq = 



X 



Rn 



i?ri 



and the coefficients depend on 

• which is just acos(l/p), 

• xi, which from (I37|) depends on p and as, 

• Xc, which from ((43]) depends on p, Xm, X and as, 



• Xm, which from (132]), only depends on p and as, 

for Ro/Rrcs < 1. 



^max/^rcs cau be Written as 



Furthermore, looking at Eq. ([24]) . for Rq ^ i?rcs, X 

as 

2 as 



X 



3K{p) 



Hp, Xr. 



1/2 



(48) 



and from Eqs. (glJ-dlT]) and (HH])-®, stiU under the 
assumption Rq <^ Rres, the terms g, h in the above 
equation depend on Rq only through as- We conclude 
that, provided that Rq <C i?rcs, equation (j47]) becomes 
frequency independent, and can in fact be written in im- 
plicit form as 



I{as,p) = 0. 



(49) 



This equation can easily be solved for as{p), in order to 
find the approximate, frequency- independent threshold 
for inversion of the Bjerknes force. The solution is pre- 
sented in the inset of Fig. Below the curve, J > 0, 
so that the Bjerknes force attracts the bubble toward 
pressure antinodes, while it becomes repulsive above. 

From as = 2a/{poRo), the inversion threshold can also 
be plotted in the (i?o,p) plane in the case of water at 



ambient pressure (cr = 0.072 N.m 



Po 



101300 Pa). 



The result is displayed in Fig. H] (thick solid line) and 
compared to the exact inversion thresholds calculated 
from numerical simulation for three driving frequencies 
20 kHz (dash-dotted line), 40 kHz (dashed line), and 
80 kHz (thin solid line). The labels on the two latter 
curves represent the value of Ro/Rrcs- It is seen that the 
above procedure yields a good estimation of the inversion 
threshold, up to Rq/R^cs — 0.1, above which it starts to 
diverge from the exact value. The reasons for this dis- 
agreement comes from the neglected Ro/Rrcs term in all 
expressions, and also from the fact that for increasing 
frequency, the bubble rebounds become more important, 
so that the bubble dynamics for x > Xc also contributes 
to expression ([M]) . Besides, a cascade of period-doubhng 
bifurcations and chaos [l^ H^l appear in some cases 
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(and are responsible for the noisy oscillations on the 80 
kHz curve) , so that the correct averaging of the Bjerknes 
force in such cases should be carried out over more than 
a single acoustic period. We did not pursue further this 
issue, since analytical predictions for these bifurcations 
are out of the scope of the present paper. 




(/iin) 

FIG. 4: Threshold of Bjerknes force inversion in the {Ro,p) 
plane, for a bubble in water in ambient conditions (u = 0.072 
N.m"\ po = 101300 Pa, fi = 10"^ Pa.s). The region / > 
corresponds to attraction by the pressure antinode, and I < 
to repulsion. The thin lines are calculated from numerical 
simulations of the RP equation. Thin solid line: / = 80 
kHz; dashed line: / = 40 kHz; dash-dotted line:/ = 20 kHz. 
The labels on the curves indicate the ratio Rq / Rrcs (triangles: 
/ = 80 kHz; filled circles: / = 40 kHz). Thick solid line: 
universal threshold calculated from approximate dynamics by 
solving (|49[) . Thick dashed line: Blake threshold. The inset 
represent the solution of (|49p in the {as,p) plane. 

Marginally, it can be seen that the inversion threshold 
in the {as,p) plane is almost linear, so that the following 
linear fit (represented by a dashed line in the inset of 
Fig. 2]) can be proposed for practical applications: 

p = 0.269 as + 1-62. (50) 

These results suggest that the inversion threshold is inde- 
pendent of frequency, and of the properties of the gas and 
liquid other than surface tension, as long as i?o/J^ros ^ 1- 
This astonishing result originates from the fact that the 
Bjerknes force mainly depends on the expansion phase of 
the bubble, which, within the approximations leading to 
Eq. (|16p , is merely governed by the driving pressure am- 
plitude and surface tension. The reasonably good agree- 
ment found in Figs. partially supports this analysis. 

In order to further investigate this issue, we first recal- 
culated the three inversion thresholds of Fig. H] (/ = 20, 
40, 80 kHz), replacing the thermal model of Ref. 32 by 
an isothermal behavior for the bubble interior. Figure [S] 



displays the results obtained (thick solid lines) and re- 
calls the thresholds calculated in Fig. 2] (thin solid lines). 
It can be seen that the thresholds slightly diverge for 
increasing Rq, but remain almost indistinguishable for 
Ro/Rrcs < 0.15 We also repeated the calculations with 
the thermal model of Ref. |32, but for argon bubbles (not 
shown), and found a negligible deviation from the air 
curves. We therefore conclude that the detailed bub- 
ble interior has a very weak influence on the expansion 
phase, at least for low enough values of Ro/Rrcs, so that 
Eq. ([50]) indeed constitutes a gas-independent law, within 
its range of validity. 



2.1 1 < r 



P 




■ 2 4 6 8 10 

FIG. 5: Same as Fig.|4l The thin solid lines are the numerical 
curves of Fig. H (/ = 20, 40, 80 kHz). The thick solid lines 
are calculated in the same conditions, except that the gas 
behavior is considered isothermal. The thin dashed line is 
the analytical threshold calculated from Eq. (|49|) . 

Another issue is the sensitivity of the results to the 
liquid viscosity. The latter has been neglected in the 
analytical approach, when approximating the RP equa- 
tion dS]) by Eq. (ITB|) . The good agreement found in Fig. |3] 
between analytical and numerical results, calculated for 
water at ambient temperature (/i = 10~^ Pa.s), suggests 
that for such low values, viscosity indeed plays a minor 
role during the bubble expansion. One should however 
check whether it is still the case for larger viscosities. 
We therefore repeated the calculation of the inversion 
threshold for viscosities 10 and 20 times larger than the 
one of water (Fig. (SJ thick dashed line and thick dash- 
dotted line) . It is clearly seen that the threshold increases 
noticeably with viscosity. Conversely, we also checked 
that the result was unaffected by decreasing the viscos- 
ity below the water's one, by computing the threshold 
for fi = 0.1/iwatcr (thick solid line). This indicates that 
viscous friction plays a non-negligible role in the bubble 
expansion for viscosities above some critical value. As al- 
ready mentioned in Ref. 1, increasing viscosity decreases 
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i?inax, and we also checked that it decreases Xm too, so 
that, strictly speaking, the Bjerknes force and its inver- 
sion threshold are viscosity-dependent. Following our re- 
sults, this influence is negligible for viscosities near or 
lower than the water's one, but for slightly larger values, 
the viscous term should be kept in the RP equation. 

2.1 1 — ^ ' ' ' ' 1 



P 




■ 2 4 6 8 10 

R() (^m) 

FIG. 6: Same as Fig.H only for / = 20 kHz, and for different 
liquid viscosities. Thin solid line: water (same as dash-dotted 
line of Fig. |4]); Thick solid line: fj, = 0.1 /^watcr; Thick dashed 
line: fi — lO/Xwater; Thick dash-dotted Une jj. — 20/iwatcr. 
The thin dashed line is the analytical threshold calculated 
from Eq. (|49)). 



V. DISCUSSION 

Important conclusions can be drawn from these re- 
sults. Figure 2] shows that the inversion thresholds for 
all frequencies (thin lines) asymptotically merge with the 
Blake threshold (thick dashed line) for small bubble radii, 
and in reasonable agreement with the analytical approx- 
imation (thick solid line). Thus, as the driving pressure 
reaches, say 1.8 bar, the range of ambient radii of iner- 
tial bubbles attracted toward the antinode is suddenly 
reduced, with an upper limit lower than 2 /xm. This ex- 
plains why a well-defined bubble-free region can be ob- 
served around t he p ressure antinode for high amplitude 
standing waves [21|. The range of attracted bubbles is 
however not void, which suggests that the zone around 
the antinode could still be filled with inertial bubbles, of 
ambient radii very close to the Blake threshold, but too 
small to be visible. As noticed in Ref. [l^, in a high am- 
plitude standing wave, the Bjerknes force acts as a sorter 
of inertial bubbles, leaving the smallest ones approaching 
or even reaching the pressure antinodes. The advantage 



of the present analysis is that it yields, through Eq. (|49l) . 
or its simpler form ([SH]) . an explicit classification of the 
bubble sizes as a function of the local acoustic pressure, 
parametrized by the ratio a/po. 

As the increasingly small bubbles approach the pres- 
sure antinode, they may coalesce or quickly grow by rec- 
tified diffusion [35*1. Increasing their size, they may again 
enter the repulsion zone in the {Ro,p) plane and move 
back again. This picture is still complicated by the poten- 
tial appearance of surface instabilities. Thus, the appar- 
ently void region observed around large pressure antin- 
odes may be in fact the locus of the complex evolution of 
very small bubbles, of sizes close to the Blake threshold. 

Finally, it is seen from the inset of Fig. |4]that decreas- 
ing as lowers the driving at which the pressure antinode 
becomes repulsive. The dimensionless parameter as can 
be varied experimentally by modifying the surface ten- 
sion a (for example adding ionic salts or surfactants), 
or by changing the static pressure Pq. The present re- 
sults suggest that, for identical bubble ambient radii, the 
Bjerknes force would become repulsive for lower drivings, 
when either decreasing a or increasing pq. This should 
have an observable effect on the size of the bubble-free 
region around the pressure antinode. However, it should 
be noted that surface tension also pla ys a crucial role 
for bubbles surface instabilities p^ . [ssl 139| , and also for 
rectified diffusion 'ssl], through the same dimensionless 
parameter as- Thus, changing as may also directly in- 
fiuence these two processes, with probable consequences 
on the bubble cloud behavior. The present result just 
demonstrates that surface tension can influence the shape 
of the bubble cloud through its direct effect on the bubble 
dynamics, and on the primary Bjerknes force. 

Figure [B] also indicates that the size of the bubble- 
free region around the pressure antinode would decrease 
noticeably when increasing viscosity slightly above the 
one of water. As mentioned in Ref. [H, this may be eas- 
ily achieved experimentally by adding glycerin in water. 
Here again, such a macroscopic effect is mediated by the 
sensitivity of the bubble dynamics to the physical prop- 
erties. To account analytically for this dependence on 
viscosity, the viscous term should be kept in the Rayleigh 
equation, which renders the approximation scheme more 
involved. A generalization of our analytical results to 
this case may be addressed in a future study. 

Finally, it is highly probable that the same effect of sur- 
face tension could be observed on the secondary Bjerknes 
force, as suggested by numerical simulations [40]. The 
extension of the present analytical method to the lat- 
ter effect is difficult, first because the expression of the 
secondary Bjerknes force also involves the bubbles veloc- 
ities, which are much more sensitive to approximations 
than the bubble radius itself, and secondly because the 
dynamics equation of the two bubbles must be coupled 
by a radiation term. 
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